Magnetic dynamics of ferromagnetic long range order in icosahedral quasicrystal

Quasicrystals lack translational symmetry and have unique lattice structures with rotational symmetry forbidden in periodic crystals. The electric state and physical property are far from complete understanding, which are the frontiers of modern matter physics. Recent discovery of the ferromagnetic long-range order in the rare-earth based icosahedral quasicrystal has made the breakthrough. Here, we first reveal the dynamical as well as static magnetic structure in the ferromagnetic long-range order in the terbium-based quasicrystal. The dynamical structure factor exhibits highly structured energy and wavenumber dependences beyond the crystalline-electric-field excitation. We find the presence of the magnetic excitation mode analog to magnon with unique hierarchical structure as well as the localized magnetic excitation with high degeneracy in the quasicrystal. Non-collinear and non-coplanar magnetic structure on the icosahedron is discovered to give rise to non-reciprocal magnetic excitation in the quasicrystal as well as non-reciprocal magnon in the periodic cubic 1/1 approximant. These findings afford illuminating insight into the magnetic dynamics in the broad range of the rare-earth-based quasicrystals and approximants.

in the mirror plane and is 80 • tilted from the pseudo 5-fold axis (see Fig. 1b) forming the ferrimagnetic state on the IC 12 . This state is uniformly distributed at the center and the corner in the unit cell of the body-center-cubic (bcc) lattice, which forms the FM-long range order in the 1/1 AC, as shown in Fig. 1c 12 .
Although the arrangement of the magnetic moments of the FM order has been identified in the real space as Fig. 1a, the magnetic structure in the reciprocal space has not been elucidated 31 . Furthermore, the dynamical magnetic structure factor has not been clarified theoretically nor experimentally. In periodic crystals, the excitation in the magnetically ordered state is known as "magnon". In the QC, the wavenumber is no longer a good quantum number and hence it is interesting whether the magnon-like excitation exists in the reciprocal space. It is also curious whether the magnetic excitation unique to the icosahedral QC exists.
So far, the lattice dynamics in the QC has been studied by inelastic X-ray-and neutron-scattering experiments 33,34 . Theoretical studies of the lattice dynamics have also been reported 33,35 . As for the magnetic dynamics, the dynamical structure factor has been calculated in the ferromagnetically aligned spins on the Fibonacci chain 36 and antiferromagnetically aligned spins on the two-dimensional octagonal tiling 37 .
However, little has been known about the magnetic dynamics in the real three-dimensional QC experimentally nor theoretically.
In this article, we for the first time clarify the magnetic dynamics of the FM long-range order in the threedimensional icosahedral QC. By analyzing the magnetic model taking into account the effect of the magnetic anisotropy arising from the CEF in the QC Au-SM-Tb, we clarify the dynamical as well as static magnetic structure. The dynamical structure factor is shown to exhibit highly structured energy and wavenumber dependences. We find that the pseudo magnon mode as well as the localized magnetic-excitation mode exists in the QC. We also find non-reciprocal magnetic excitation in the QC as well as the non-reciprocal magnon in the 1/1 AC. These findings provide insight into the understanding of the recently discovered FM long-range order in the QC Au 65 Ga 20 Tb 15 13 .

Results
Lattice structure of quasicrystal. Let us start with the lattice structure of the QC. Although the detailed lattice structure of the QC Au 65 Ga 20 Tb 15 has not been solved experimentally, the Tb lattice is considered to form the Cd 5.7 Yb-type QC 3 . Figure 1a shows the main structure of the QC where the Tb-12 cluster i.e., the IC is located at 30 vertices of the icosidodecahedron. In the Cd 5.7 Yb-type QC, there exist a few other ICs as well as Tb sites located between ICs. In this study, as a first step of analysis, we consider the Tb sites on ICs shown in Fig. 1a with the total lattice number is N = 12 × 30 = 360 . In this study, we employ the real Tb configuration of the IC determined in the 1/1 AC Au 70 Si 17 Tb 13 12 and locate them at the 30 vertices of the τ 3 -times enlarged icosidodecahedron in the Tsai-type cluster of Au 70 Si 17 Tb 13 . Here, τ = (1 + √ 5)/2 is the golden mean.  where J ij is the exchange interaction between the "spin" S i on the ith Tb site and S j on the jth Tb site. Here S i expresses the total angular momentum J . Since the Hund's rule tells us that the ground multiplet of the Tb 3+ ion with 4f 8 configuration is given by J = 6 , we set S = 6 . We consider the nearest neighbor (N.N.) interaction J 1 and the next N.N. (N.N.N.) interaction J 2 not only for the intra IC but also for the inter IC as discussed in Ref. 31 . The second term in Eq. (1) represents the magnetic anisotropy arising from the CEF, where e 3 is the unit vector set to be along the direction of the ordered moment at each Tb site (see Fig. 1d) and D is the parameter of anisotropy.
In the strong anisotropy limit of the model (1), it has been confirmed that the FM long-range order of the ferrimagnetic state shown in Fig. 1a is stabilized in the wide parameter region of J 2 /J 1 for 64.1 • ≤ θ ≤ 80 • of the ground-state phase diagram (see Fig. S3a in Supplementary Information) 31 . In this study, we analyze the FM ground state of the ferrimagnetic state with θ = 80 • (see Fig. 1b) in the QC lattice shown in Fig. 1a under the open boundary condition.

Magnetic excitation in quasicrystal.
Since the ferrimagnetic state is the non-collinear alignment of "spins", it is convenient to introduce the local coordinate at each Tb site where the ê 3 axis is taken as the ordered "spin" direction as shown in Fig. 1d (see "Methods" section). Then, by applying the Holstein-Primakoff transformation 38 to H, the "spin" operators are transformed to the boson operators as is the raising (lowering) "spin" operator and a † i (a i ) is a creation (annifiration) operator of the boson at the ith Tb site. We retain the quadratic terms of the boson operators since the higher order terms are considered to be irrelevant at least for the ground state. By diagonalizing H, we obtain the energy spectrum ω i for J 1 = −1 and J 2 = −1 as shown in Fig. 2a. Fig. 2b) reflecting the FM long-range order of the ferrimagnetic state for θ = 80 • (see Fig. 1a,b) as the stable ground state. We find that there appear several gaps, among which the energy gaps are visible as the step structures at i = 300 and i = 330 (see red and blue arrows in Fig. 2a respectively). Interestingly, we find that ω i from i = 234 to 263 are degenerate, i.e., degeneracy is 30, irrespective of D. It turns out that these degenerated states give rise to unique magnetic excitation in the dynamical structure factor, which will be discussed later. Below we show the results for D = 10 as a representative case for the QC Au-SM-Tb. For D = 10 , the lowest energy is ω 360 /(|J 1 |S) = 9.15 ≡ � which is defined as the gap and the largest energy is ω 1 /(|J 1 |S) = 19.57 . The degenerated energy mentioned above is ω i /(|J 1 |S) = 15.00 for i = 234, . . . 263.
Static structure factor of magnetism. To clarify the magnetic structure in the reciprocal space, we calculate the magnetic structure factor The result of F s (q) in the q z −q x plane for q y = 0 is shown in Fig. 3a. We find that a sharp peak appears at q = 0 . At finite q , sharp peaks also appear with smaller intensities. The top view is shown in Fig. 3b, where the bright spots are located along the pseudo 5-fold axis indicated by the dashed line named the d e * 2 line with an arrow. As shown in Fig. 3e, d * e i (i = 1, . . . , 6) is the primitive vector of the six-dimensional reciprocal lattice space as www.nature.com/scientificreports/ the physical (external) space components 39 . The slope of the d e * 2 line is 1.736 reflecting the real configuration of the Tb sites in the IC 12 employed in this study and it is noted that the slope is to be τ in the case of the regular IC 39 . In Fig. 3b, the bright spots are also located along the dashed line with the negative slope − 1.736, which is named the d e * 3 line. We also plot the top views of F s (q) in the q y −q z plane for q x = 0 as Fig. 3c and in the q x −q y plane for q z = 0 as Fig. 3d, where d e * i lines (i = 4, 6 and i = 1, 5) are drawn. From these results, we confirmed that the largest peak is located at q = 0 . This indicates that the uniform long-range order of the ferrimagnetic state on the IC (see Fig. 1a) is realized.
Dynamical structure factor of magnetism. The dynamical magnetic structure factor is defined as Here, |GS� is the ground state with E 0 being the ground-state energy and we set η = 10 −6 . We have calculated S αα (q, ω) (α = x, y, z) and below we show the results of S yy (q, ω) as a representative for the dynamical structure factor (see Supplementary Information for the results of S xx (q, ω) and S zz (q, ω)).
The elastic component of the dynamical structure factor S yy (q, ω = 0) is shown for q along the d e * 2 line with q z ∈ [0, 2] Å −1 in Fig. 4a. The largest peak appears at q = 0 and several peak structures also appear for finite q , whose values are typically O(10 8 ) . As shown in Fig. 2a, in the magnetic excitation, the gap is caused by the uniaxial anisotropy D arising from the CEF. To grasp the ω dependence of the dynamical structure factor, we plot S yy (q = 0, ω) in Fig. 4b. For ω/(|J 1 |S) ≥ � , the spiky peak structure appears where the peak values are typically O(10 4 ∼ 10 5 ).
Next, inelastic part of the dynamical structure factor S yy (q, ω) above the CEF excitation gap is shown in Fig. 4c for q along the d e * 2 line with q z ∈ [0, 2] Å −1 . The result shows the highly structured energy and wavenumber dependences. The largest peak appears at (q, ω/(|J 1 |S)) = (0, �) . The sharp peak structures appear in the lower energy regions above the CEF excitation gap: � ≤ ω/(|J 1 |S) < 12.7 . Interestingly, it is remarkable that the successive mountain-like high intensity structure with a large periodicity q ∼ 0.6 Å −1 appears at ω/(|J 1 |S) = 15.0 . From the relation of wavenumber and wavelength �q = 2π/ , the wavelength is estimated to be ∼ 0.67a Å. This roughly corresponds to the diameter of the IC d = 0.72a as shown in Fig. 1d. At ω/(|J 1 |S) = 15.0 , as |q| increases, the intensity decreases with this periodicity. This mode is completely localized, which appears at the flat branch, i.e., non-dispersive constant ω with 30 degeneracy, reflecting the degenerated energy ω i from i = 234   To analyze the wavenumber dependence of the magnetic excitation in more detail, we plot the top view of S yy (q, ω) in Fig. 4d. We find that successive magnon-like mode i.e. sinusoidal-like shape with periodicity about q ∼ 0.15 Å −1 appears around 12.0 < ω/(|J 1 |S) < 12.7 . These sinusoidal-like modes further form the larger sinusoidal-like modes with q ∼ 0.3 Å −1 around 12.0 < ω/(|J 1 |S) < 13.3 . These series of excitations form unique hierarchical structure. The recursive structure i.e., self-similar structure was also reported in Fibonacci chain [36] and in two-dimensional octagonal tiling 37 . The similar self-similar structure can also be seen for the lower excitation energies � < ω/(|J 1 |S) < 10.5 . We refer to these modes as pseudo-magnon modes hereafter. The pseudo magnon mode propagating along the 5-fold direction shown in Fig. 4d is inherent in the QC.
At ω/(|J 1 |S) = 15.0 , the flat branch appears, which is the localized mode as remarkably seen in Fig. 4c. Just below this branch, bright intensities exist and above the localized mode, the broad continuum-like structures appear for 15.0 < ω/(|J 1 |S) < 19.6.
As shown in Fig. 4c,d, the high intensity appears in the low-energy region for ω/(|J 1 |S) ≤ 15.0 and in the high-energy region, the low intensity appears. This is characteristic of the magnetic excitation in the FM order of the ferrimagnetic state (see Fig. 1a,b) irrespective of the strength of the anisotropy. Namely, the similar feature was confirmed to appear for 1.8 ≤ D ≤ 30.0.

Comparison with 1/1 approximant crystal. To get insight into the emergence of the pseudo magnon mode
in the QC, we calculate the magnon dispersion in the 1/1 AC. Namely, the model (1) is applied to the 1/1 AC (see Fig. 1c). On the atomic positions, here we employ the lattice structure of the 1/1 AC Au 70 Si 17 Tb 13 12 with the lattice constant being a = 14.726 Å. For the same parameter used in the QC, i.e., J 1 = −1 , J 2 = −1 , and D = 10 with S = 6 , we confirmed that the FM long-range order of the ferrimagnetic state (Fig. 1b) is realized, as shown in Fig. 1c. Then, by applying the linear spin-wave theory 38 to the model (1), we calculate the excitation energy of the magnon ω q in the reciprocal space. The result for q along the symmetry line in the bcc Brillouin zone is plotted as solid lines in Fig. 5a, where the number of the unit cell is taken as 64 3 . Similarly to the QC, owing to the uniaxial anisotropy D = 10 , the energy gap opens for 0 < ω q /(|J 1 |S) < 10 in the magnon excitation. For 10 ≤ ω q /(|J 1 |S) ≤ 25.3 , the magnon excitation forms the dispersive energy bands. Here, for the comparison with the QC, the energy range of � ≤ ω q /(|J 1 |S) ≤ 19.57 is shown in Fig. 5a (for the result of the whole energy range, see Supplementary Information). It is noted that the length of the reciprocal lattice between the Ŵ and H points is 2π/a = 0.43Å −1 .
For comparison, we calculate S yy (q, ω) in the QC for the same q as that shown in Fig. 5a. The result for � ≤ ω q /(|J 1 |S) ≤ 19.57 is plotted as the top view in Fig. 5b. We see that the dispersive modes with the broad intensity start from the Ŵ point at ω/(|J 1 |S) = � and also at ω/(|J 1 |S) = 12.5 , which seem to correspond to the dispersive magnon bands starting from the Ŵ point at ω q /(|J 1 |S) = 10 and 13.3 in Fig. 5a, respectively. These two magnon bands have the concave dispersions along the P-Ŵ -N line in Fig. 5a, which seem to correspond to   Fig. 5b. The quadratic dispersion ω = cq 2 starting from the Ŵ point at ω = � is evaluated by using the high-intensity data (red data in Fig. 5b) in the QC. The magnon velocity v = ∂ω/(∂q) = 2cq is about 4-times smaller than that estimated from the lowest-excitation-energy data in the vicinity of the Ŵ point shown in Fig. 5a in the 1/1 AC. Furthermore, the convex magnon band at the H point around ω q /(|J 1 |S) = 17.5 in Fig. 5a seems to correspond to the mountain-like broad intensities at the H point around ω/(|J 1 |S) = 17 in Fig. 5b. These results indicate that even in the QC, the magnetic excitation mode like magnon, i.e., the pseudo magnon mode, exists with broad width in their intensities. Fig. 5a, we also plotted ω −q as dashed lines. Notable is that there exists a remarkable difference between the solid line and the dashed line indicating ω q = ω −q . This implies that the non-reciprocal magnon appears in the present FM long-range order of the ferrimagnetic state. This is, to our best knowledge, the first discovery of the nonreciprocal magnon in the AC. Note that in the bcc lattice, the spatial inversion symmetry exists. The emergence of the non-reciprocal magnon is ascribed to the non-collinear "spin" configuration shown in Fig. 1b. To check this point, we calculated the magnon dispersion for the collinear "spin" configuration such as the FM long-range order where all the "spins" are aligned to the same direction in the 1/1 AC. In this case, ω q = ω −q holds for all q , which indicates that the reciprocal magnon is realized (see Supplementary Information).

Non-reciprocal magnon in 1/1 approximant crystal. In the calculation of the magnon dispersion in
Non-reciprocal magnetic excitation in quasicrystal. Next, we calculate the dynamical structure factor S yy (q, ω) and also S yy (−q, ω) in the QC for q along the symmetry lines in the bcc Brillouin zone of the 1/1 AC, as shown in Fig. 6a,b, respectively. By comparing both results, we find that there exist several differences between q and −q , i.e., S yy (q, ω) = S yy (−q, ω) . For instance, near the Ŵ point around ω/(|J 1 |S) ≈ 17 , a clear difference exists. Indeed, the non-reciprocal magnon dispersion, i.e., deviation from the dashed line from the solid line, is shown to appear near the Ŵ point around ω q /(|J 1 |S) ≈ 17 in Fig. 5a. Furthermore, around the H  www.nature.com/scientificreports/ point for 17 ≤ ω/(|J 1 |S) ≤ 19 , remarkable differences exist as seen in Fig. 6a,b. Actually, the ω q = ω −q feature is prominent at the H point around ω q /(|J 1 |S) ≈ 17 in Fig. 5a. This is, to our best knowledge, the first discovery of the non-reciprocal magnetic excitation in the FM long-range order in the QC. This finding suggests that the non-reciprocal magnetic excitation appears generally in the icosahedral QC. This is because the alignment of the "spins" at the 12 vertices of the IC (see Fig. 1b) following the uniaxial anisotropy induced by the CEF at each rare-earth site inevitably becomes the non-collinear type.

Discussion.
We have clarified the static and dynamical structure factor in the FM long-range order of the ferrimagnetic state in the icosahedral QC. The static structure factor exhibits the largest peak at q = 0 , which indicates the uniform long-range order of the ferrimagnetic state. By calculating the dynamical structure factor S(q, ω) , the elastic part for ω = 0 is shown to have the largest peak at q = 0 . The uniaxial anisotropy arising from the CEF causes the energy gap in the magnetic excitation so that the non-elastic part exhibits the intensity, whose order is 4 decades smaller than that of the elastic part, beyond the CEF gap. The inelastic S(q, ω) shows the highly structured energy and wavenumber dependences. In spite of no translational invariance in the QC, we have found the pseudo magnon mode with hierarchical structure propagating along the 5-fold direction inherent in the QC. We have also found the completely localized magnetic-excitation mode with periodicity characterized by the wavelength corresponding to the size of the IC in the wavenumber space as the degenerated modes. This implies that this mode reflects the magnetic excitation on each icosahedron because of the cancellation of the dispersions of the excitation arising from symmetric interaction J 1 = J 2 (see section II in Supplementary  Information). The non-collinear and non-coplanar "spin" structure on the IC is shown to give rise to the nonreciprocal magnetic excitation in the QC as well as the non-reciprocal magnon in the AC. Our model and the results are expected to be relevant to the broad range of the rare-earth based icosahedral QCs and ACs with strong magnetic anisotropy. Hence, our findings as well as the method for the analysis developed in this study open a new research field of the magnetic dynamics in the QCs and ACs. It is noted that consideration of more detailed lattice structures of the Cd 5.7 Yb-type QC and the examination of the effect of the system size i.e., boundary condition on magnetic dynamics are left for future studies.

Methods
Quasicrystal and approximant crystal. The rare-earth based QC and AC consists of the Tsai-type cluster with nested shell structures of polyhedrons. The rare-earth atom is located at the 12 vertices of the IC. The AC retains the periodicity as well as the common local atomic configuration to the QC. There exists a series of the ACs such as 1/1 AC, 2/1 AC, 3/2 AC, . . . , where F n in the F n+1 /F n AC is the Fibonacci number (F n+2 = F n+1 + F n , F 1 = F 2 = 1) . In the 1/1 AC, the IC composed of the rare-earth atoms is distributed at the center and corner of the bcc lattice. As n increases, the size of the unit cell of the F n+1 /F n AC expands and for the n → ∞ limit (lim n→∞ F n+1 /F n = τ ) the size of the unit cell becomes infinite, which corresponds to the QC.

Theory of magnetic excitation in quasicrystal.
To calculate the magnetic excitation from the FM long-range order in the QC, we transform the spin operators in the model (1) into the boson operators. Since the ferrimagnetic state is a noncolinear magnetic state, it is convenient to introduce the local coordinate at each Tb site 40 . The unit vectors in the global xyz coordinate r 1 =x , r 2 =ŷ , and r 3 =ẑ are expressed by the local orthogonal coordinate with the unit vector ê i 3 , whose direction is indicated by the polar angles (θ i , φ i ) , as (see Fig. 1d). Here, R i is the rotation matrix defined as Then, the first term in Eq. (1) is expressed as By using S i ·ê i 1 = (S + i + S − i )/2 and S i ·ê i 2 = (S + i − S − i )/(2i) where S + i and S − i are raising and lowering "spin" operators, respectively, we apply the Holstein-Primakoff transformation 38 to H. Namely, "spin" operators are expressed by the boson operators as S + i = √ 2S − n i a i , S − i = a † i √ 2S − n i and S i ·ê i 3 = S − n i with n i ≡ a † i a i . We retain the quadratic terms with respect to a † i and a i , which are considered to be at least valid for the ground state. In the noncollinear magnetic state as the hedgehog, anomalous terms such as a † i a † j and a i a j appear. The resultant H is expressed as where χ † = (a † 1 , a † 2 , . . . , a † N ) and is the 2N × 2N matrix. By performing the para unitary transformation (4) r α = R i αβê i β (5) R i = cos θ i cos φ i − sin φ i sin θ i cos φ i cos θ i sin φ i cos φ i sin θ i sin φ i − sin θ i 0 cos θ i .

Data availability
All the data supporting the findings are available from the corresponding author upon reasonable request.